Benchmarking Models 3 - Informed Priors
March 18, 2022
Introduction
Performs a Bayesian calculation with uninformed priors to estimate
model 2 of
\[ Y_{icb2} = \alpha_{b} + \delta^{T} T_{c} + \delta^{GK} T^{GK}_c + \beta X_{icb1} + \rho Y_{icb1} +\gamma_{1} \tau_{c} +\epsilon_{icb2} \]
Warnings (if applicable)
Setup
the needful
Load the relevant packages
Set global core settings, if applicable. Use global setting only if necessary, best to use as an inline setting to avoid over allocating system resources
The models - table 3
These models build out table 3 which presents the results for the primary and secondary outcomes utilizing informed priors with 1. normal distributions and 2. Cauchy distributions.We follow the advice of Rachael Meager and incorporate disagreement in the literature and use 6 standard deviations and a mean = 0
Per Capita Monthly Consumption - Primary Outcome
Load the Monthly per capita consumption model variables
Normal Distribution
per capita consumption: This is the basic benchmarking model utilizing the a global normal distribution with mean = 0 and 6 standard deviations.
Model Summery
summary(per_cap_consumption_normal_bayesmodel)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: consumption_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + consumption_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + (1 | block) + (1 | vid)
## Data: per_cap_consumption_data (Number of observations: 1750)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.27 0.06 0.16 0.40 1.00 925 1840
##
## ~vid (Number of levels: 248)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.38 0.03 0.32 0.44 1.00 1822 2416
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 8.51 0.18 8.16 8.87 1.00 4268
## cost_deviation 0.00 0.00 0.00 0.00 1.00 2799
## treat_any 0.11 0.08 -0.05 0.28 1.00 2023
## treat_GK -0.13 0.08 -0.30 0.03 1.00 2313
## consumption_asinh_R1 0.18 0.01 0.16 0.21 1.00 6550
## Lhh_wealth_asinh 0.02 0.01 0.01 0.03 1.00 8404
## Lvill_eligible_ratio 0.16 0.33 -0.48 0.80 1.00 1889
## Tail_ESS
## Intercept 3113
## cost_deviation 2941
## treat_any 2617
## treat_GK 3076
## consumption_asinh_R1 3162
## Lhh_wealth_asinh 2465
## Lvill_eligible_ratio 2542
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 1.16 0.01 1.14 1.19 1.00 6660 2748
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(per_cap_consumption_normal_bayesmodel)## prior class coef group resp dpar nlpar
## normal(0,6) b
## normal(0,6) b consumption_asinh_R1
## normal(0,6) b cost_deviation
## normal(0,6) b Lhh_wealth_asinh
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## student_t(3, 10.7, 2.5) Intercept
## student_t(3, 0, 2.5) sd
## student_t(3, 0, 2.5) sd block
## student_t(3, 0, 2.5) sd Intercept block
## student_t(3, 0, 2.5) sd vid
## student_t(3, 0, 2.5) sd Intercept vid
## student_t(3, 0, 2.5) sigma
## bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(per_cap_consumption_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept uninformative
## 2 b_cost_deviation uninformative
## 3 b_treat_any uninformative
## 4 b_treat_GK uninformative
## 5 b_consumption_asinh_R1 uninformative
## 6 b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio uninformativeDiagnostics
# trace diagnostic plot
mcmc_trace(per_cap_consumption_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_consumption_asinh_R1", "b_Lhh_wealth_asinh",
"b_Lvill_eligible_ratio", "sd_block__Intercept",
"sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\per_cap_consumption_norm_trace.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plot
mcmc_dens(per_cap_consumption_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_consumption_asinh_R1", "b_Lhh_wealth_asinh",
"b_Lvill_eligible_ratio", "sd_block__Intercept",
"sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\per_cap_consumption_norm_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(per_cap_consumption_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_consumption_asinh_R1", "b_Lhh_wealth_asinh",
"b_Lvill_eligible_ratio", "sd_block__Intercept",
"sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\per_cap_consumption_norm_overlay.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens(per_cap_consumption_normal_bayesmodel,pars = c("b_treat_GK"))
mcmc_dens_overlay(per_cap_consumption_normal_bayesmodel,pars = c("b_treat_GK"))
#acf (auto-correlation) diagnostic plot
mcmc_acf(per_cap_consumption_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_consumption_asinh_R1", "b_Lhh_wealth_asinh",
"b_Lvill_eligible_ratio", "sd_block__Intercept",
"sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\per_cap_consumption_norm_acf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(per_cap_consumption_normal_bayesmodel, nsamples = 100)## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.
pp_check(per_cap_consumption_normal_bayesmodel, nsamples = 10, type = 'error_scatter_avg', alpha = .1)## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.
Dietary Diversity - Primary Outcome
Load the Dietary Diversity variables
Dietary Diversity: This is the basic bechmarking model utilizing the default, uninformed priors
dietary_diversity_normal_bayesmodel <-
brm(formula = dietarydiversity | weights(samp_wgt) ~
cost_deviation + treat_any + treat_GK +
dietarydiversity_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lsavingsstock_asinh3 +
Lconsumpti_x_Ldietarydi + Lconsumpti_x_Lproductiv + Ldietarydi_x_Lassetscon +
(1 | vid) + (1 | block),
data = dietary_diversity_data,
family = gaussian("identity"),
set_prior("normal(0,6)", class = "b"),
seed = 1272022,
warmup = 1000,
iter = 2000,
thin = 1,
control = list(adapt_delta = .95, max_treedepth = 10),
#backend = "cmdstanr",
cores = 4, #overrides default 1 core
#threads = 3,need to get cmdstanr package working here
save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
file = "informed_prior_outcomes\\dietary_diversity_normal_bayes")# tidy_dietary_diversity_normal_bayesmodel <- tidy(dietary_diversity_normal_bayesmodel)
# #view(tidy_dietary_diversity_normal_bayesmodel)
# write_csv(tidy_dietary_diversity_normal_bayesmodel, "informed_prior_outcomes\\dietary_diversity_normal_bayes.csv")Model Summary
summary(dietary_diversity_normal_bayesmodel)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: dietarydiversity | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + dietarydiversity_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lsavingsstock_asinh3 + Lconsumpti_x_Ldietarydi + Lconsumpti_x_Lproductiv + Ldietarydi_x_Lassetscon + (1 | vid) + (1 | block)
## Data: dietary_diversity_data (Number of observations: 1751)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.30 0.08 0.16 0.48 1.01 999 1490
##
## ~vid (Number of levels: 248)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.55 0.04 0.47 0.64 1.00 1524 1977
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 3.19 0.22 2.75 3.62 1.00 3506
## cost_deviation 0.00 0.00 0.00 0.00 1.00 2913
## treat_any 0.25 0.12 0.02 0.48 1.00 2215
## treat_GK -0.03 0.11 -0.26 0.19 1.00 1985
## dietarydiversity_R1 -0.12 0.07 -0.27 0.02 1.00 3123
## Lhh_wealth_asinh 0.01 0.01 -0.00 0.03 1.00 5196
## Lvill_eligible_ratio -0.57 0.44 -1.46 0.30 1.00 1857
## Lsavingsstock_asinh3 0.00 0.00 0.00 0.00 1.00 4079
## Lconsumpti_x_Ldietarydi 0.01 0.01 0.00 0.03 1.00 3153
## Lconsumpti_x_Lproductiv 0.01 0.00 0.00 0.01 1.00 4034
## Ldietarydi_x_Lassetscon 0.01 0.00 0.01 0.01 1.00 4947
## Tail_ESS
## Intercept 2970
## cost_deviation 2827
## treat_any 2894
## treat_GK 2504
## dietarydiversity_R1 2880
## Lhh_wealth_asinh 2539
## Lvill_eligible_ratio 2270
## Lsavingsstock_asinh3 3187
## Lconsumpti_x_Ldietarydi 2679
## Lconsumpti_x_Lproductiv 3445
## Ldietarydi_x_Lassetscon 3143
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 1.58 0.02 1.55 1.62 1.00 5780 3160
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(dietary_diversity_normal_bayesmodel)## prior class coef group resp dpar nlpar
## normal(0,6) b
## normal(0,6) b cost_deviation
## normal(0,6) b dietarydiversity_R1
## normal(0,6) b Lconsumpti_x_Ldietarydi
## normal(0,6) b Lconsumpti_x_Lproductiv
## normal(0,6) b Ldietarydi_x_Lassetscon
## normal(0,6) b Lhh_wealth_asinh
## normal(0,6) b Lsavingsstock_asinh3
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## student_t(3, 5, 2.5) Intercept
## student_t(3, 0, 2.5) sd
## student_t(3, 0, 2.5) sd block
## student_t(3, 0, 2.5) sd Intercept block
## student_t(3, 0, 2.5) sd vid
## student_t(3, 0, 2.5) sd Intercept vid
## student_t(3, 0, 2.5) sigma
## bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(dietary_diversity_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept uninformative
## 2 b_cost_deviation uninformative
## 3 b_treat_any uninformative
## 4 b_treat_GK uninformative
## 5 b_dietarydiversity_R1 uninformative
## 6 b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio uninformative
## 8 b_Lsavingsstock_asinh3 uninformative
## 9 b_Lconsumpti_x_Ldietarydi uninformative
## 10 b_Lconsumpti_x_Lproductiv uninformative
## 11 b_Ldietarydi_x_Lassetscon uninformativeDiagnostics
# trace diagnostic plot
mcmc_trace(dietary_diversity_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_dietarydiversity_R1", "b_Lhh_wealth_asinh",
"b_Lvill_eligible_ratio", "b_Lsavingsstock_asinh3",
"b_Lconsumpti_x_Ldietarydi", "b_Lconsumpti_x_Lproductiv",
"b_Ldietarydi_x_Lassetscon", "sd_block__Intercept",
"sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\dietary_div_normal_trace.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plot
mcmc_dens(dietary_diversity_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_dietarydiversity_R1", "b_Lhh_wealth_asinh",
"b_Lvill_eligible_ratio", "b_Lsavingsstock_asinh3",
"b_Lconsumpti_x_Ldietarydi", "b_Lconsumpti_x_Lproductiv",
"b_Ldietarydi_x_Lassetscon", "sd_block__Intercept",
"sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\dietary_div_normal_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(dietary_diversity_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_dietarydiversity_R1", "b_Lhh_wealth_asinh",
"b_Lvill_eligible_ratio", "b_Lsavingsstock_asinh3",
"b_Lconsumpti_x_Ldietarydi", "b_Lconsumpti_x_Lproductiv",
"b_Ldietarydi_x_Lassetscon", "sd_block__Intercept",
"sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\dietary_div_normal_overlay.png", plot = last_plot(), width = 12, height = 5)
#acf (auto-correlation) diagnostic plot
mcmc_acf(dietary_diversity_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_dietarydiversity_R1", "b_Lhh_wealth_asinh",
"b_Lvill_eligible_ratio", "b_Lsavingsstock_asinh3",
"b_Lconsumpti_x_Ldietarydi", "b_Lconsumpti_x_Lproductiv",
"b_Ldietarydi_x_Lassetscon", "sd_block__Intercept",
"sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\dietary_div_normal_acf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(dietary_diversity_normal_bayesmodel, nsamples = 100)## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.
pp_check(dietary_diversity_normal_bayesmodel, nsamples = 10, type = 'error_scatter_avg', alpha = .1)## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.
Total Household Wealth - Primary Outcome
Load the total household wealth model variables
total household wealth: This is the basic
benchmarking model utilizing brm() default, uninformed
priors. Removed Lhh_wealth_asinh to account for collinearity issues.
hh_wealth_normal_bayesmodel <-
brm(formula = wealth_asinh | weights(samp_wgt) ~
cost_deviation + treat_any + treat_GK +
wealth_asinh_R1 + Lvill_eligible_ratio + Lowndwelling +
(1 | block) + (1 | vid),
data = hh_wealth_data,
family = gaussian("identity"),
set_prior("normal(0,6)", class = "b"),
seed = 1272022,
warmup = 1000,
iter = 2000,
thin = 1,
control = list(adapt_delta = .99, max_treedepth = 10),
#backend = "cmdstanr",
cores = 4, #overrides default 1 core
#threads = 3,need to get cmdstanr package working here
save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
file = "informed_prior_outcomes\\hh_wealth_normal_bayes")Model Summery
summary(hh_wealth_normal_bayesmodel)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: wealth_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + wealth_asinh_R1 + Lvill_eligible_ratio + Lowndwelling + (1 | block) + (1 | vid)
## Data: hh_wealth_data (Number of observations: 1751)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.52 0.15 0.24 0.84 1.00 838 1059
##
## ~vid (Number of levels: 248)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 1.03 0.09 0.85 1.22 1.00 1546 2483
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 7.81 0.32 7.18 8.42 1.00 3221
## cost_deviation 0.00 0.00 -0.00 0.00 1.00 3080
## treat_any 0.06 0.23 -0.40 0.52 1.00 2187
## treat_GK 0.01 0.23 -0.43 0.47 1.00 2311
## wealth_asinh_R1 0.18 0.02 0.15 0.21 1.00 5499
## Lvill_eligible_ratio -0.10 0.84 -1.79 1.49 1.00 1882
## Lowndwelling 3.46 0.20 3.07 3.86 1.00 5247
## Tail_ESS
## Intercept 3327
## cost_deviation 3530
## treat_any 2541
## treat_GK 2728
## wealth_asinh_R1 3585
## Lvill_eligible_ratio 2751
## Lowndwelling 3514
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 3.43 0.04 3.35 3.51 1.00 6564 2972
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(hh_wealth_normal_bayesmodel)## prior class coef group resp dpar nlpar
## normal(0,6) b
## normal(0,6) b cost_deviation
## normal(0,6) b Lowndwelling
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## normal(0,6) b wealth_asinh_R1
## student_t(3, 13.9, 2.5) Intercept
## student_t(3, 0, 2.5) sd
## student_t(3, 0, 2.5) sd block
## student_t(3, 0, 2.5) sd Intercept block
## student_t(3, 0, 2.5) sd vid
## student_t(3, 0, 2.5) sd Intercept vid
## student_t(3, 0, 2.5) sigma
## bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(hh_wealth_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept uninformative
## 2 b_cost_deviation uninformative
## 3 b_treat_any uninformative
## 4 b_treat_GK uninformative
## 5 b_wealth_asinh_R1 uninformative
## 6 b_Lvill_eligible_ratio informative
## 7 b_Lowndwelling uninformativeDiagnostics
# trace diagnostic plot
mcmc_trace(hh_wealth_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_wealth_asinh_R1",
"b_Lvill_eligible_ratio", "b_Lowndwelling",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\hh_wealth_normal_trace.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plot
mcmc_dens(hh_wealth_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_wealth_asinh_R1",
"b_Lvill_eligible_ratio", "b_Lowndwelling",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\hh_wealth_normal_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(hh_wealth_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_wealth_asinh_R1",
"b_Lvill_eligible_ratio", "b_Lowndwelling",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\hh_wealth_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)
#acf (auto-correlation) diagnostic plot
mcmc_acf(hh_wealth_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any",
"b_treat_GK", "b_wealth_asinh_R1",
"b_Lvill_eligible_ratio", "b_Lowndwelling",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\hh_wealth_normal_acf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(hh_wealth_normal_bayesmodel, nsamples = 100)## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.
pp_check(hh_wealth_normal_bayesmodel, nsamples = 10, type = 'error_scatter_avg', alpha = .1)## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.
Borrowing Stock - Secondary Outcome
Load the borrowing stock model variables
borrowing stock: This is the basic benchmarking
model utilizing brm() global weakly informed priors with
mean = 0 and 6 standard deviations.
borrowing_stock_normal_bayesmodel <-
brm(formula = borrowingstock_asinh | weights(samp_wgt) ~
cost_deviation + treat_any + treat_GK +
borrowingstock_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio +
(1 | block) + (1 | vid),
data = borrowing_stock_data,
family = gaussian("identity"),
set_prior("normal(0,6)", class = "b"),
seed = 1272022,
warmup = 1000,
iter = 2000,
thin = 1,
control = list(adapt_delta = .95, max_treedepth = 10),
#backend = "cmdstanr",
cores = 4, #overrides default 1 core
#threads = 3,need to get cmdstanr package working here
save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
file = "informed_prior_outcomes\\borrowing_stock_normal_bayes")Model Summery
summary(borrowing_stock_normal_bayesmodel)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: borrowingstock_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + borrowingstock_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + (1 | block) + (1 | vid)
## Data: borrowing_stock_data (Number of observations: 1751)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.93 0.20 0.58 1.37 1.00 1115 1692
##
## ~vid (Number of levels: 248)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 1.19 0.11 0.99 1.41 1.00 1551 2673
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 6.13 0.43 5.28 6.97 1.00 2114
## cost_deviation 0.00 0.00 -0.00 0.00 1.00 3212
## treat_any -0.73 0.28 -1.28 -0.19 1.00 2114
## treat_GK 0.67 0.28 0.13 1.23 1.00 1993
## borrowingstock_asinh_R1 0.23 0.02 0.20 0.26 1.00 4909
## Lhh_wealth_asinh -0.00 0.02 -0.04 0.03 1.00 4630
## Lvill_eligible_ratio -0.50 1.08 -2.56 1.69 1.00 1775
## Tail_ESS
## Intercept 2566
## cost_deviation 3023
## treat_any 2525
## treat_GK 2440
## borrowingstock_asinh_R1 3375
## Lhh_wealth_asinh 3116
## Lvill_eligible_ratio 2416
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 4.48 0.05 4.38 4.59 1.00 6601 2872
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(borrowing_stock_normal_bayesmodel)## prior class coef group resp dpar nlpar
## normal(0,6) b
## normal(0,6) b borrowingstock_asinh_R1
## normal(0,6) b cost_deviation
## normal(0,6) b Lhh_wealth_asinh
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## student_t(3, 9.2, 2.9) Intercept
## student_t(3, 0, 2.9) sd
## student_t(3, 0, 2.9) sd block
## student_t(3, 0, 2.9) sd Intercept block
## student_t(3, 0, 2.9) sd vid
## student_t(3, 0, 2.9) sd Intercept vid
## student_t(3, 0, 2.9) sigma
## bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(borrowing_stock_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept uninformative
## 2 b_cost_deviation uninformative
## 3 b_treat_any uninformative
## 4 b_treat_GK uninformative
## 5 b_borrowingstock_asinh_R1 uninformative
## 6 b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio informativeDiagnostics
# trace diagnostic plot
mcmc_trace(borrowing_stock_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_borrowingstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\borrowing_stock_normal_trace.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plots
mcmc_dens(borrowing_stock_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_borrowingstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\borrowing_stock_normal_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(borrowing_stock_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_borrowingstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\borrowing_stock_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)
#acf (auto-correlation) diagnostic plot
mcmc_acf(borrowing_stock_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_borrowingstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\borrowing_stock_normal_acf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(borrowing_stock_normal_bayesmodel, nsamples = 100)## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.
pp_check(borrowing_stock_normal_bayesmodel, nsamples = 10, type = 'error_scatter_avg', alpha = .1)## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.
Savings Stock - Secondary Outcome
Load the savings stock model variables
savings stock: This is the basic bechmarking model
utilzing brm() default uninformed priors
savings_stock_normal_bayesmodel <-
brm(formula = savingsstock_asinh | weights(samp_wgt) ~
cost_deviation + treat_any + treat_GK +
savingsstock_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio +
Lconsumpti_x_Lproductiv + Lconsumpti_x_Lassetscon +
(1 | block) + (1 | vid),
data = savings_stock_data,
family = gaussian("identity"),
set_prior("normal(0,6)", class = "b"),
seed = 1272022,
warmup = 1000,
iter = 2000,
thin = 1,
control = list(adapt_delta = .95, max_treedepth = 10),
#backend = "cmdstanr",
cores = 4, #overrides default 1 core
#threads = 3,need to get cmdstanr package working here
save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
file = "informed_prior_outcomes\\savings_stock_normal_bayes")Model Summery
summary(savings_stock_normal_bayesmodel)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: savingsstock_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + savingsstock_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lconsumpti_x_Lproductiv + Lconsumpti_x_Lassetscon + (1 | block) + (1 | vid)
## Data: savings_stock_data (Number of observations: 1751)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.70 0.19 0.36 1.12 1.00 840 1205
##
## ~vid (Number of levels: 248)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 1.40 0.11 1.18 1.63 1.01 1074 1983
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 1.78 0.46 0.88 2.68 1.00 2347
## cost_deviation 0.00 0.00 -0.00 0.00 1.00 2022
## treat_any -0.13 0.30 -0.73 0.48 1.00 1409
## treat_GK 1.31 0.30 0.73 1.92 1.00 1287
## savingsstock_asinh_R1 0.17 0.01 0.14 0.20 1.00 4002
## Lhh_wealth_asinh -0.05 0.02 -0.09 -0.02 1.00 4623
## Lvill_eligible_ratio 1.54 1.17 -0.68 3.89 1.00 1232
## Lconsumpti_x_Lproductiv 0.02 0.00 0.02 0.03 1.00 5177
## Lconsumpti_x_Lassetscon 0.01 0.00 0.01 0.01 1.00 5780
## Tail_ESS
## Intercept 3166
## cost_deviation 2700
## treat_any 1854
## treat_GK 2084
## savingsstock_asinh_R1 2628
## Lhh_wealth_asinh 3062
## Lvill_eligible_ratio 1905
## Lconsumpti_x_Lproductiv 3323
## Lconsumpti_x_Lassetscon 3550
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 4.15 0.05 4.06 4.25 1.00 5102 3097
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(savings_stock_normal_bayesmodel)## prior class coef group resp dpar nlpar
## normal(0,6) b
## normal(0,6) b cost_deviation
## normal(0,6) b Lconsumpti_x_Lassetscon
## normal(0,6) b Lconsumpti_x_Lproductiv
## normal(0,6) b Lhh_wealth_asinh
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b savingsstock_asinh_R1
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## student_t(3, 8.7, 3.1) Intercept
## student_t(3, 0, 3.1) sd
## student_t(3, 0, 3.1) sd block
## student_t(3, 0, 3.1) sd Intercept block
## student_t(3, 0, 3.1) sd vid
## student_t(3, 0, 3.1) sd Intercept vid
## student_t(3, 0, 3.1) sigma
## bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(savings_stock_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept uninformative
## 2 b_cost_deviation uninformative
## 3 b_treat_any uninformative
## 4 b_treat_GK uninformative
## 5 b_savingsstock_asinh_R1 uninformative
## 6 b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio informative
## 8 b_Lconsumpti_x_Lproductiv uninformative
## 9 b_Lconsumpti_x_Lassetscon uninformativeDiagnostics
# trace diagnostic plot
mcmc_trace(savings_stock_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_savingsstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lconsumpti_x_Lproductiv", "b_Lconsumpti_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\savings_stock_normal_trace.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plots
mcmc_dens(savings_stock_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_savingsstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lconsumpti_x_Lproductiv", "b_Lconsumpti_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\savings_stock_normal_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(savings_stock_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_savingsstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lconsumpti_x_Lproductiv", "b_Lconsumpti_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\savings_stock_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)
#acf (auto-correlation) diagnostic plot
mcmc_acf(savings_stock_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_savingsstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lconsumpti_x_Lproductiv", "b_Lconsumpti_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\savings_stock_normal_acf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(savings_stock_normal_bayesmodel, nsamples = 100)## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.
pp_check(savings_stock_normal_bayesmodel, nsamples = 10, type = 'error_scatter_avg', alpha = .1)## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.
Health Knowledge - Secondary Outcome
Load the health knowledge model variables
health knowledge: This is the basic bechmarking
model utilzing brm() default, uninformed priors
health_knowledge_normal_bayesmodel <-
brm(formula = health_knowledge | weights(samp_wgt) ~
cost_deviation + treat_any + treat_GK +
health_knowledge_R1 +
Lhh_wealth_asinh + Lvill_eligible_ratio +
(1 | block) + (1 | vid),
data = health_knowledge_data,
family = gaussian("identity"),
set_prior("normal(0,6)", class = "b"),
seed = 1272022,
warmup = 1000,
iter = 2000,
thin = 1,
control = list(adapt_delta = .95, max_treedepth = 10),
#backend = "cmdstanr",
cores = 4, #overrides default 1 core
#threads = 3,need to get cmdstanr package working here
save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
file = "informed_prior_outcomes\\health_knowledge_normal_bayes")Model Summery
summary(health_knowledge_normal_bayesmodel)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: health_knowledge | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + health_knowledge_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + (1 | block) + (1 | vid)
## Data: health_knowledge_data (Number of observations: 1751)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.34 0.22 0.02 0.81 1.00 454 1098
##
## ~vid (Number of levels: 248)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 1.64 0.12 1.42 1.88 1.00 1474 2879
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 1.37 0.40 0.60 2.17 1.00 2668
## cost_deviation 0.00 0.00 -0.00 0.00 1.00 2071
## treat_any 0.10 0.33 -0.54 0.76 1.00 1447
## treat_GK -0.08 0.34 -0.76 0.60 1.00 1364
## health_knowledge_R1 0.06 0.02 0.03 0.09 1.00 4833
## Lhh_wealth_asinh 0.11 0.02 0.07 0.15 1.00 5603
## Lvill_eligible_ratio 0.61 1.15 -1.66 2.89 1.00 1349
## Tail_ESS
## Intercept 2726
## cost_deviation 2465
## treat_any 2351
## treat_GK 2339
## health_knowledge_R1 2690
## Lhh_wealth_asinh 2820
## Lvill_eligible_ratio 1861
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 4.35 0.05 4.25 4.45 1.00 4968 2742
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(health_knowledge_normal_bayesmodel)## prior class coef group resp dpar nlpar
## normal(0,6) b
## normal(0,6) b cost_deviation
## normal(0,6) b health_knowledge_R1
## normal(0,6) b Lhh_wealth_asinh
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## student_t(3, 3.2, 3.8) Intercept
## student_t(3, 0, 3.8) sd
## student_t(3, 0, 3.8) sd block
## student_t(3, 0, 3.8) sd Intercept block
## student_t(3, 0, 3.8) sd vid
## student_t(3, 0, 3.8) sd Intercept vid
## student_t(3, 0, 3.8) sigma
## bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(health_knowledge_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept informative
## 2 b_cost_deviation uninformative
## 3 b_treat_any uninformative
## 4 b_treat_GK uninformative
## 5 b_health_knowledge_R1 uninformative
## 6 b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio informativeDiagnostics
# trace diagnostic plot
mcmc_trace(health_knowledge_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_health_knowledge_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\health_knowledge_normal_trace.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plots
mcmc_dens(health_knowledge_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_health_knowledge_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\health_knowledge_normal_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(health_knowledge_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_health_knowledge_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\health_knowledge_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)
#acf (auto-correlation) diagnostic plot
mcmc_acf(health_knowledge_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_health_knowledge_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\health_knowledge_normal_acf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(health_knowledge_normal_bayesmodel, ndraws = 100)
pp_check(health_knowledge_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)
sanitation practices - Secondary Outcome
Load the sanitation practices model variables
sanitation practices: This is the basic benchmarking model utilizing the default, uninformed priors
sanitation_practices_normal_bayesmodel <-
brm(formula = sanitation_practices | weights(samp_wgt) ~
cost_deviation + treat_any + treat_GK +
sanitation_practices_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio +
Lproductiv_x_Lassetscon +
(1 | block) + (1 | vid),
data = sanitation_practices_data,
family = gaussian("identity"),
set_prior("normal(0,6)", class = "b"),
seed = 1272022,
warmup = 1000,
iter = 2000,
thin = 1,
control = list(adapt_delta = .95, max_treedepth = 10),
#backend = "cmdstanr",
cores = 4, #overrides default 1 core
#threads = 3,need to get cmdstanr package working here
save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
file = "informed_prior_outcomes\\sanitation_practices_normal_bayes")Model Summery
summary(sanitation_practices_normal_bayesmodel)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: sanitation_practices | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + sanitation_practices_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lproductiv_x_Lassetscon + (1 | block) + (1 | vid)
## Data: sanitation_practices_data (Number of observations: 1751)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.18 0.12 0.01 0.44 1.00 695 1527
##
## ~vid (Number of levels: 248)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.97 0.07 0.84 1.10 1.00 1975 3150
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept -2.01 0.24 -2.49 -1.54 1.00 2838
## cost_deviation 0.00 0.00 -0.00 0.00 1.00 2594
## treat_any 0.21 0.20 -0.18 0.61 1.00 1824
## treat_GK -0.31 0.20 -0.70 0.09 1.00 1809
## sanitation_practices_R1 0.11 0.02 0.07 0.14 1.00 6618
## Lhh_wealth_asinh 0.04 0.01 0.01 0.06 1.00 8011
## Lvill_eligible_ratio 0.03 0.69 -1.32 1.41 1.00 1698
## Lproductiv_x_Lassetscon 0.01 0.00 0.01 0.01 1.00 5045
## Tail_ESS
## Intercept 3352
## cost_deviation 2614
## treat_any 2550
## treat_GK 2471
## sanitation_practices_R1 3076
## Lhh_wealth_asinh 3119
## Lvill_eligible_ratio 2384
## Lproductiv_x_Lassetscon 3458
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 2.63 0.03 2.57 2.69 1.00 7096 2932
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(sanitation_practices_normal_bayesmodel)## prior class coef group resp dpar nlpar
## normal(0,6) b
## normal(0,6) b cost_deviation
## normal(0,6) b Lhh_wealth_asinh
## normal(0,6) b Lproductiv_x_Lassetscon
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b sanitation_practices_R1
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## student_t(3, 0.1, 3) Intercept
## student_t(3, 0, 3) sd
## student_t(3, 0, 3) sd block
## student_t(3, 0, 3) sd Intercept block
## student_t(3, 0, 3) sd vid
## student_t(3, 0, 3) sd Intercept vid
## student_t(3, 0, 3) sigma
## bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(sanitation_practices_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept informative
## 2 b_cost_deviation uninformative
## 3 b_treat_any uninformative
## 4 b_treat_GK uninformative
## 5 b_sanitation_practices_R1 uninformative
## 6 b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio informative
## 8 b_Lproductiv_x_Lassetscon uninformativeDiagnostics
# trace diagnostic plot
mcmc_trace(sanitation_practices_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_sanitation_practices_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lproductiv_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\sanitation_practices_normal_trace.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plots
mcmc_dens(sanitation_practices_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_sanitation_practices_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lproductiv_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\sanitation_practices_normal_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(sanitation_practices_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_sanitation_practices_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lproductiv_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\sanitation_practices_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)
#acf (auto-correlation) diagnostic plot
mcmc_acf(sanitation_practices_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_sanitation_practices_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lproductiv_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\sanitation_practices_normal_acf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(sanitation_practices_normal_bayesmodel, ndraws = 100)
pp_check(sanitation_practices_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)
productive assets - Secondary Outcome
Load the productive assets model variables
productive assets: This is the basic benchmarking model utilizing the default, uninformed priors
productive_assets_normal_bayesmodel <-
brm(formula = productiveassets_asinh | weights(samp_wgt) ~
cost_deviation + treat_any + treat_GK +
productiveassets_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio +
Lconsumpti_x_Lassetscon +
(1 | block) + (1 | vid),
data = productive_assets_data,
family = gaussian("identity"),
set_prior("normal(0,6)", class = "b"),
seed = 1272022,
warmup = 1000,
iter = 2000,
thin = 1,
control = list(adapt_delta = .95, max_treedepth = 10),
#backend = "cmdstanr",
cores = 4, #overrides default 1 core
#threads = 3,need to get cmdstanr package working here
save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
file = "informed_prior_outcomes\\productive_assets_normal_bayes")Model Summery
summary(productive_assets_normal_bayesmodel)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: productiveassets_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + productiveassets_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lconsumpti_x_Lassetscon + (1 | block) + (1 | vid)
## Data: productive_assets_data (Number of observations: 1751)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.34 0.07 0.22 0.50 1.00 1918 2377
##
## ~vid (Number of levels: 248)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.43 0.04 0.36 0.51 1.00 1525 2705
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 6.89 0.21 6.48 7.31 1.00 4710
## cost_deviation 0.00 0.00 0.00 0.00 1.00 4761
## treat_any 0.37 0.10 0.18 0.57 1.00 3640
## treat_GK -0.27 0.10 -0.46 -0.07 1.00 3640
## productiveassets_asinh_R1 0.33 0.02 0.30 0.37 1.00 9267
## Lhh_wealth_asinh -0.00 0.01 -0.01 0.01 1.00 8694
## Lvill_eligible_ratio 0.12 0.40 -0.67 0.88 1.00 3211
## Lconsumpti_x_Lassetscon 0.01 0.00 0.00 0.01 1.00 4417
## Tail_ESS
## Intercept 3456
## cost_deviation 3691
## treat_any 3388
## treat_GK 3483
## productiveassets_asinh_R1 3190
## Lhh_wealth_asinh 3238
## Lvill_eligible_ratio 3084
## Lconsumpti_x_Lassetscon 3279
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 1.50 0.02 1.47 1.53 1.00 7656 3088
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(productive_assets_normal_bayesmodel)## prior class coef group resp dpar
## normal(0,6) b
## normal(0,6) b cost_deviation
## normal(0,6) b Lconsumpti_x_Lassetscon
## normal(0,6) b Lhh_wealth_asinh
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b productiveassets_asinh_R1
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## student_t(3, 11.4, 2.5) Intercept
## student_t(3, 0, 2.5) sd
## student_t(3, 0, 2.5) sd block
## student_t(3, 0, 2.5) sd Intercept block
## student_t(3, 0, 2.5) sd vid
## student_t(3, 0, 2.5) sd Intercept vid
## student_t(3, 0, 2.5) sigma
## nlpar bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(productive_assets_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept uninformative
## 2 b_cost_deviation uninformative
## 3 b_treat_any uninformative
## 4 b_treat_GK uninformative
## 5 b_productiveassets_asinh_R1 uninformative
## 6 b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio uninformative
## 8 b_Lconsumpti_x_Lassetscon uninformativeDiagnostics
# trace diagnostic plot
mcmc_trace(productive_assets_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_productiveassets_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lconsumpti_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\productive_assets_normal_trace.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plots
mcmc_dens(productive_assets_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_productiveassets_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lconsumpti_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\productive_assets_normal_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(productive_assets_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_productiveassets_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lconsumpti_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\productive_assets_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)
#acf (auto-correlation) diagnostic plot
mcmc_acf(productive_assets_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_productiveassets_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lconsumpti_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\productive_assets_normal_cf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(productive_assets_normal_bayesmodel, ndraws = 100)
pp_check(productive_assets_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)
Consumption Assets - Secondary Outcome
Load the consumption assets model variables
consumption assets: This is the basic benchmarking model utilizing a global normal prior with mean = and 6 standard deviations. This model has a number of errors: Warning:
- There were 4 chains where the estimated Bayesian Fraction of Missing Information was low. See https://mc-stan.org/misc/warnings.html#bfmi-low
- Warning: Examine the pairs() plot to diagnose sampling problems
- Warning: The largest R-hat is 4.48, indicating chains have not mixed. Running the chains for more iterations may help. See https://mc-stan.org/misc/warnings.html#r-hat
- Warning: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable.
- Warning: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable. https://mc-stan.org/misc/warnings.html#tail-ess
consumption_assets_normal_bayesmodel <-
brm(formula = assetsconsumption_asinh | weights(samp_wgt) ~
cost_deviation + treat_any + treat_GK +
assetsconsumption_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb + Ldurablesexpenditure +
Ldietarydi_x_Lassetscon + Lproductiv_x_Lassetscon +
(1 | block) + (1 | vid),
data = consumption_assets_data,
family = gaussian("identity"),
set_prior("normal(0,6)", class = "b"),
seed = 1272022,
warmup = 1000,
iter = 2000,
thin = 1,
control = list(adapt_delta = .95, max_treedepth = 10),
#backend = "cmdstanr",
cores = 4, #overrides default 1 core
#threads = 3,need to get cmdstanr package working here
save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
file = "informed_prior_outcomes\\consumption_assets_normal_bayes")Model Summery
summary(consumption_assets_normal_bayesmodel)## Warning: Parts of the model have not converged (some Rhats are > 1.05). Be
## careful when analysing the results! We recommend running more iterations and/or
## setting stronger priors.## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: assetsconsumption_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + assetsconsumption_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb + Ldurablesexpenditure + Ldietarydi_x_Lassetscon + Lproductiv_x_Lassetscon + (1 | block) + (1 | vid)
## Data: consumption_assets_data (Number of observations: 1751)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.89 1.09 0.00 2.72 3.49 4 14
##
## ~vid (Number of levels: 248)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.82 0.55 0.33 1.71 4.43 4 11
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept -0.35 5.35 -15.70 4.06 3.42 4
## cost_deviation 0.00 0.00 0.00 0.01 2.12 5
## treat_any 0.54 1.22 -1.47 1.88 3.51 4
## treat_GK 0.93 1.07 -1.23 1.79 3.56 4
## assetsconsumption_asinh_R1 0.34 0.35 -0.63 0.82 2.46 5
## Lhh_wealth_asinh -0.17 0.23 -0.57 -0.01 1.60 7
## Lvill_eligible_ratio 0.70 0.50 -0.15 1.36 3.28 5
## Lroomsnumb -0.16 0.81 -1.49 0.95 2.93 5
## Ldurablesexpenditure -0.00 0.00 -0.00 0.00 1.72 6
## Ldietarydi_x_Lassetscon 0.18 0.29 0.01 0.70 1.74 6
## Lproductiv_x_Lassetscon -0.01 0.07 -0.18 0.08 2.72 5
## Tail_ESS
## Intercept 11
## cost_deviation 33
## treat_any 11
## treat_GK 13
## assetsconsumption_asinh_R1 11
## Lhh_wealth_asinh 15
## Lvill_eligible_ratio 15
## Lroomsnumb 11
## Ldurablesexpenditure 11
## Ldietarydi_x_Lassetscon 11
## Lproductiv_x_Lassetscon 11
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 7.29 6.27 3.33 19.92 2.53 6 13
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(consumption_assets_normal_bayesmodel)## prior class coef group resp dpar
## normal(0,6) b
## normal(0,6) b assetsconsumption_asinh_R1
## normal(0,6) b cost_deviation
## normal(0,6) b Ldietarydi_x_Lassetscon
## normal(0,6) b Ldurablesexpenditure
## normal(0,6) b Lhh_wealth_asinh
## normal(0,6) b Lproductiv_x_Lassetscon
## normal(0,6) b Lroomsnumb
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## student_t(3, 9.8, 2.5) Intercept
## student_t(3, 0, 2.5) sd
## student_t(3, 0, 2.5) sd block
## student_t(3, 0, 2.5) sd Intercept block
## student_t(3, 0, 2.5) sd vid
## student_t(3, 0, 2.5) sd Intercept vid
## student_t(3, 0, 2.5) sigma
## nlpar bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(consumption_assets_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept informative
## 2 b_cost_deviation uninformative
## 3 b_treat_any informative
## 4 b_treat_GK informative
## 5 b_assetsconsumption_asinh_R1 uninformative
## 6 b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio uninformative
## 8 b_Lroomsnumb informative
## 9 b_Ldurablesexpenditure uninformative
## 10 b_Ldietarydi_x_Lassetscon uninformative
## 11 b_Lproductiv_x_Lassetscon uninformativeDiagnostics
# trace diagnostic plot
mcmc_trace(consumption_assets_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_assetsconsumption_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lroomsnumb", "b_Ldurablesexpenditure",
"b_Ldietarydi_x_Lassetscon", "b_Lproductiv_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\consumption_assets_normal_race.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plots
mcmc_dens(consumption_assets_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_assetsconsumption_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lroomsnumb", "b_Ldurablesexpenditure",
"b_Ldietarydi_x_Lassetscon", "b_Lproductiv_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\consumption_assets_normal_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(consumption_assets_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_assetsconsumption_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lroomsnumb", "b_Ldurablesexpenditure",
"b_Ldietarydi_x_Lassetscon", "b_Lproductiv_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\consumption_assets_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)
#acf (auto-correlation) diagnostic plot
mcmc_acf(consumption_assets_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_assetsconsumption_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lroomsnumb", "b_Ldurablesexpenditure",
"b_Ldietarydi_x_Lassetscon", "b_Lproductiv_x_Lassetscon",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\consumption_assets_normal_acf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(consumption_assets_normal_bayesmodel, ndraws = 100)
pp_check(consumption_assets_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)
House Value - Secondary Outcome
Load the house value model variables
house value: This is the basic benchmarking model utilizing a prior normal distribution with mean = 0 and 6 standard deviations. This model had a number of errors and the posterior estimates are unreliable.
- Warning: There were 2023 divergent transitions after warmup. See https://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup to find out why this is a problem and how to eliminate them.
- Warning: There were 1977 transitions after warmup that exceeded the maximum treedepth. Increase max_treedepth above 10. See https://mc-stan.org/misc/warnings.html#maximum-treedepth-exceeded
- Warning: There were 4 chains where the estimated Bayesian Fraction of Missing Information was low. See https://mc-stan.org/misc/warnings.html#bfmi-low
- Warning: Examine the pairs() plot to diagnose sampling problems
- Warning: The largest R-hat is 4.46, indicating chains have not mixed. Running the chains for more iterations may help. See https://mc-stan.org/misc/warnings.html#r-hat
- Warning: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable. Running the chains for more iterations may help. See https://mc-stan.org/misc/warnings.html#bulk-ess
- Warning: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable. Running the chains for more iterations may help. See https://mc-stan.org/misc/warnings.html#tail-ess
dwelling_cost_normal_bayesmodel <-
brm(formula = selfcostdwell_asinh | weights(samp_wgt) ~
cost_deviation + treat_any + treat_GK +
selfcostdwell_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb + Ldurablesexpenditure +
(1 | block) + (1 | vid),
data = dwelling_cost_data,
family = gaussian("identity"),
set_prior("normal(0,6)", class = "b"),
seed = 1272022,
warmup = 1000,
iter = 2000,
thin = 1,
control = list(adapt_delta = .95, max_treedepth = 10),
#backend = "cmdstanr",
cores = 4, #overrides default 1 core
#threads = 3,need to get cmdstanr package working here
save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
file = "informed_prior_outcomes\\dwelling_cost_normal_bayes")Model Summery
summary(dwelling_cost_normal_bayesmodel)## Warning: Parts of the model have not converged (some Rhats are > 1.05). Be
## careful when analysing the results! We recommend running more iterations and/or
## setting stronger priors.## Warning: There were 2023 divergent transitions after warmup. Increasing
## adapt_delta above 0.95 may help. See http://mc-stan.org/misc/
## warnings.html#divergent-transitions-after-warmup## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: selfcostdwell_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + selfcostdwell_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb + Ldurablesexpenditure + (1 | block) + (1 | vid)
## Data: dwelling_cost_data (Number of observations: 1654)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.83 0.66 0.20 1.93 3.32 4 11
##
## ~vid (Number of levels: 247)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.29 0.15 0.16 0.70 4.15 4 11
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 12.56 8.41 -1.15 26.73 3.35 4
## cost_deviation 0.00 0.00 -0.00 0.00 2.59 5
## treat_any 0.13 1.02 -1.45 1.51 2.44 5
## treat_GK 1.51 0.27 0.94 1.79 3.17 4
## selfcostdwell_asinh_R1 -0.20 0.86 -1.33 0.80 2.56 5
## Lhh_wealth_asinh 0.07 0.06 -0.01 0.19 2.51 5
## Lvill_eligible_ratio 0.69 0.53 -0.19 1.28 3.42 4
## Lroomsnumb 0.06 0.98 -1.43 1.31 3.07 5
## Ldurablesexpenditure 0.00 0.00 0.00 0.00 1.85 6
## Tail_ESS
## Intercept 11
## cost_deviation 30
## treat_any 22
## treat_GK 11
## selfcostdwell_asinh_R1 31
## Lhh_wealth_asinh 22
## Lvill_eligible_ratio 11
## Lroomsnumb 11
## Ldurablesexpenditure 21
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 3.06 1.95 0.89 8.63 2.91 5 11
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(dwelling_cost_normal_bayesmodel)## prior class coef group resp dpar nlpar
## normal(0,6) b
## normal(0,6) b cost_deviation
## normal(0,6) b Ldurablesexpenditure
## normal(0,6) b Lhh_wealth_asinh
## normal(0,6) b Lroomsnumb
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b selfcostdwell_asinh_R1
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## student_t(3, 13.8, 2.5) Intercept
## student_t(3, 0, 2.5) sd
## student_t(3, 0, 2.5) sd block
## student_t(3, 0, 2.5) sd Intercept block
## student_t(3, 0, 2.5) sd vid
## student_t(3, 0, 2.5) sd Intercept vid
## student_t(3, 0, 2.5) sigma
## bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(dwelling_cost_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept informative
## 2 b_cost_deviation uninformative
## 3 b_treat_any informative
## 4 b_treat_GK uninformative
## 5 b_selfcostdwell_asinh_R1 informative
## 6 b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio uninformative
## 8 b_Lroomsnumb informative
## 9 b_Ldurablesexpenditure uninformativeDiagnostics
# trace diagnostic plot
mcmc_trace(dwelling_cost_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_selfcostdwell_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lroomsnumb", "b_Ldurablesexpenditure",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\dwelling_cost_normal_trace.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plots
mcmc_dens(dwelling_cost_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_selfcostdwell_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lroomsnumb", "b_Ldurablesexpenditure",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\dwelling_cost_normal_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(dwelling_cost_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_selfcostdwell_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lroomsnumb", "b_Ldurablesexpenditure",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\dwelling_cost_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)
#acf (auto-correlation) diagnostic plot
mcmc_acf(dwelling_cost_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_selfcostdwell_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
"b_Lroomsnumb", "b_Ldurablesexpenditure",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\dwelling_cost_normal_acf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(dwelling_cost_normal_bayesmodel, ndraws = 100)
pp_check(dwelling_cost_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)
house quality - Secondary Outcome
Load the house quality model variables
house quality: This is the basic benchmarking model utilizing the default, uninformed priors
housing_quality_normal_bayesmodel <-
brm(formula = housing_quality | weights(samp_wgt) ~
cost_deviation + treat_any + treat_GK +
housing_quality_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb +
(1 | block) + (1 | vid),
data = housing_quality_data,
family = gaussian("identity"),
set_prior("normal(0,6)", class = "b"),
seed = 1272022,
warmup = 1000,
iter = 2000,
thin = 1,
control = list(adapt_delta = .95, max_treedepth = 10),
#backend = "cmdstanr",
cores = 4, #overrides default 1 core
#threads = 3,need to get cmdstanr package working here
save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
file = "informed_prior_outcomes\\housing_quality_normal_bayes")Model Summery
summary(housing_quality_normal_bayesmodel)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: housing_quality | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + housing_quality_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb + (1 | block) + (1 | vid)
## Data: housing_quality_data (Number of observations: 1751)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~block (Number of levels: 22)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.11 0.07 0.01 0.27 1.01 635 1399
##
## ~vid (Number of levels: 248)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.59 0.05 0.50 0.68 1.00 2097 2836
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept -2.62 0.21 -3.02 -2.21 1.00 4181
## cost_deviation 0.00 0.00 -0.00 0.00 1.00 2716
## treat_any -0.12 0.13 -0.38 0.14 1.00 2334
## treat_GK -0.08 0.13 -0.34 0.18 1.00 2093
## housing_quality_R1 0.01 0.02 -0.03 0.05 1.00 5163
## Lhh_wealth_asinh 0.04 0.01 0.02 0.05 1.00 5704
## Lvill_eligible_ratio 0.14 0.43 -0.69 0.99 1.00 1752
## Lroomsnumb 0.52 0.04 0.45 0.59 1.00 4986
## Tail_ESS
## Intercept 3366
## cost_deviation 2751
## treat_any 2695
## treat_GK 2609
## housing_quality_R1 3668
## Lhh_wealth_asinh 2826
## Lvill_eligible_ratio 2517
## Lroomsnumb 2988
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 1.92 0.02 1.88 1.96 1.00 6705 3269
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Prior summery - how informative are priors
prior_summary(housing_quality_normal_bayesmodel)## prior class coef group resp dpar nlpar
## normal(0,6) b
## normal(0,6) b cost_deviation
## normal(0,6) b housing_quality_R1
## normal(0,6) b Lhh_wealth_asinh
## normal(0,6) b Lroomsnumb
## normal(0,6) b Lvill_eligible_ratio
## normal(0,6) b treat_any
## normal(0,6) b treat_GK
## student_t(3, -0.1, 2.5) Intercept
## student_t(3, 0, 2.5) sd
## student_t(3, 0, 2.5) sd block
## student_t(3, 0, 2.5) sd Intercept block
## student_t(3, 0, 2.5) sd vid
## student_t(3, 0, 2.5) sd Intercept vid
## student_t(3, 0, 2.5) sigma
## bound source
## user
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## default
## default
## (vectorized)
## (vectorized)
## (vectorized)
## (vectorized)
## defaultcheck_prior(housing_quality_normal_bayesmodel)## Parameter Prior_Quality
## 1 b_Intercept informative
## 2 b_cost_deviation uninformative
## 3 b_treat_any uninformative
## 4 b_treat_GK uninformative
## 5 b_housing_quality_R1 uninformative
## 6 b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio uninformative
## 8 b_Lroomsnumb uninformativeDiagnostics
# trace diagnostic plot
mcmc_trace(housing_quality_normal_bayesmodel, n_warmup = 0,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_housing_quality_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", "b_Lroomsnumb",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\housing_quality_normal_trace.png", plot = last_plot(), width = 12, height = 5)
#density diagnostic plots
mcmc_dens(housing_quality_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_housing_quality_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", "b_Lroomsnumb",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\housing_quality_normal_dens.png", plot = last_plot(), width = 12, height = 5)
mcmc_dens_overlay(housing_quality_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_housing_quality_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", "b_Lroomsnumb",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\housing_quality_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)
#acf (auto-correlation) diagnostic plot
mcmc_acf(housing_quality_normal_bayesmodel,
pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
"b_housing_quality_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", "b_Lroomsnumb",
"sd_block__Intercept", "sd_vid__Intercept", "sigma"))
ggsave("table_3_diagnostics\\housing_quality_normal_acf.png", plot = last_plot(), width = 12, height = 5)posterior predictive checks
pp_check(housing_quality_normal_bayesmodel, ndraws = 100)
pp_check(housing_quality_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)